Effects of grazing intensity on diversity and composition of rhizosphere and non‐rhizosphere microbial communities in a desert grassland

Abstract Overgrazing‐induced grassland degradation has become a serious ecological problem worldwide. The diversity and composition of soil microbial communities are sensitive to grazing disturbances. However, our understanding is limited with respect to the effects of grazing intensity on bacterial and fungal communities, especially in plant rhizosphere. Using a long‐term grazing experiment, we evaluated the diversity and composition of microbial communities in both rhizosphere and non‐rhizosphere soils under three grazing intensities (light, moderate, and heavy grazing) in a desert grassland and examined the relative roles of grazing‐induced changes in some abiotic and biotic factors in affecting the diversity and composition of microbial communities. Our results showed that soil bacteria differed greatly in diversity and composition between rhizosphere and non‐rhizosphere zones, and so did soil fungi. Moderate and heavy grazing significantly reduced the rhizosphere bacterial diversity. Grazing intensity substantially altered the bacterial composition and the fungal composition in both zones but with different mechanisms. While root nitrogen and soil nitrogen played an important role in shaping the rhizosphere bacterial composition, soil‐available phosphorus greatly affected the non‐rhizosphere bacterial composition and the fungal composition in both soils. This study provides direct experimental evidence that the diversity and composition of microbial communities were severely altered by heavy grazing on a desert grassland. Thus, to restore the grazing‐induced, degraded grasslands, we should pay more attention to the conservation of soil microbes in addition to vegetation recovery.


| INTRODUC TI ON
Grasslands cover about 40% of the earth's land surface and are mainly used for livestock grazing (Suttie et al., 2005). While light and/or moderate grazing can improve grassland diversity and productivity, heavy grazing generally causes severe grassland degradation (Kemp et al., 2013;McNaughton, 1979). By far, overgrazing-induced grassland degradation has become a major ecological problem worldwide (Dlamini et al., 2016;Eldridge et al., 2017;van de Koppel et al., 1997). This is especially the case in the desert steppe ecosystems on the Mongolian Plateau, where the vegetation cover is ca. 25% and the annual biomass production is ca. 130 g m −2 (Kemp et al., 2013;Wei et al., 2013).
Overgrazing not only reduces vegetation cover, decreases plant production, and alters vegetation composition (Wei et al., 2013) but also has profound impacts on the diversity and composition of soil microbial communities (Hu et al., 2017;Lamb et al., 2011;Wang, Ding, et al., 2021;. Some studies in desert steppe ecosystems showed that grazing could significantly reduce the diversity of microbial communities in bulk soil and in rhizosphere soil (Wang, Ding, et al., 2021;Zhang et al., 2019). In contrast, another study in the same ecosystem showed that grazing had little effect on bacterial diversity and fungal diversity in rhizosphere soil, nor in bulk soil . These inconsistent results suggest that the impacts of grazing on microbial communities are likely to be mediated by grazing intensity and grazing history. Moreover, a previous study indicated that grazing could alter the composition of microbial communities in rhizosphere soil rather than in non-rhizosphere soil , highlighting the contrasting responses of microbial communities between rhizosphere and non-rhizosphere soils. In fact, soil microbial communities can differ greatly in diversity and composition between the two zones because host plants usually exert a select effect on microbial taxa in their rhizosphere soil (Hartmann et al., 2009;Nan et al., 2020). Considering that grazing has different impacts on different plant species and that different plant species prefer different microbial taxa, the impacts of grazing on rhizosphere microbial communities are expected to be stronger than on non-rhizosphere microbial communities. However, how grazing intensity affects the soil bacterial and fungal communities in rhizosphere and non-rhizosphere zones, especially over a decade, remains poorly understood.
Grazing can influence the soil microbe through several mechanisms. First, herbivory on plants can have dual effects on soil microbial communities. On the one hand, removal of aboveground biomass decreases the input of litter into the soil, which reduces the supply of carbon resource to the soil microbe, and thus causes a reduction in soil microbial richness (Dwivedi et al., 2017). On the other hand, removal of aboveground parts can also stimulate the allocation of organic nutrients to belowground organs (Mueller et al., 2017), resulting in enrichment of root exudates (Hamilton & Frank, 2001), which in turn affects the structure of microbial community . Second, manure from livestock can affect soil pH values (Wang, Ding, et al., 2021;Yang, Niu, et al., 2018;Yang, Zhu, et al., 2018) and improve the available soil nitrogen (Han et al., 2008). As a result, the relative abundance of copiotrophic taxa was enhanced, while that of oligotrophic taxa was reduced (Nie et al., 2018). Finally, livestock trampling decreases soil moisture due to increased soil compactness (Kobayashi et al., 1997), which in turn affects soil microbial properties . These studies provide important insights into how grazing-induced changes in biotic and abiotic factors can influence microbial communities; however, the relative roles of biotic and abiotic factors in affecting microbial responses to grazing intensity in rhizosphere and non-rhizosphere soils are poorly evaluated.
Given that desert steppe is more sensitive to grazing disturbance than other steppes (e.g., typical steppe or meadow steppe) due to its relatively lower precipitation, lower vegetation cover, and lower biomass productivity (Wei et al., 2013), we addressed these research gaps by a 15-year grazing intensity experiment in a desert steppe. We focused on three specific questions: (i) How do microbial communities differ in diversity and composition between rhizosphere and non-rhizosphere? (ii) How do the bacterial and fungal communities respond to grazing intensities? (iii) What are the relative roles of abiotic and biotic factors in affecting the responses of bacterial and fungal communities to grazing intensity? Overall, this study contributes to our understanding of the impacts of grazing intensity on the diversity and composition of microbial communities in rhizosphere and non-rhizosphere soils, as well as the mechanisms driving these changes.

| Study site
This study was carried out in a desert grassland on the Inner Mongolia Plateau (41°46′43.6″ N, 111°53′41.7″ E; elevation 1456 m). This area has a temperate continental grassland climate, with an annual mean temperature and precipitation of 3.6°C and 214 mm, respectively (Wang et al., 2014). About 80% of precipitation occurs during the growing season (May to September). The vegetation at the experiment site was dominated by Stipa breviflora and Cleistogenes songorica. The soil is a Kastanozem soil (FAO soil classification) with a sandy loam texture. The Siziwang Banner, Inner Mongolia Autonomous Region of China, is traditionally used as pastures for free grazing and 90% of them have been degraded due to overgrazing (Kemp et al., 2013).

| Experimental design
In 2004, a long-term grazing intensity experiment was established. We had four grazing intensities: control (zero sheep-unit ha −1 month −1 , G0), light (0.15 sheep-unit ha −1 month −1 , G0.15), moderate (0.30 sheep-unit ha −1 month −1 , G0.30), and heavy grazing (0.45 sheep-unit ha −1 month −1 , G0.45; Zhang et al., 2021). G0.15, G0.30, and G0.45 were approximately corresponding to 18%, 40%, and 60% of annual biomass consumed by animals, respectively. Each treatment had three replicates. Under a complete randomized block design, 12 plots, each with an area of 4.4 ha, were separated by fences. Sheep of similar size were selected in advance and the grazing period lasted from June 1 to November 30 every year. The plots were grazed from June 1 to November 30 every year and from 6 AM to 6 PM every day. At night, the sheep were placed in an enclosure.

| Field sampling and laboratory analysis
We randomly selected 30 bunches of C. songorica in each plot in late August 2018, corresponding to the peak biomass time. Each bunch of plants together with soil attached to their roots were dug out to a depth of 15 cm (Cao et al., 2019;Yang, Niu, et al., 2018;Yang, Zhu, et al., 2018), with approximately 90% of the roots harvested. Nonrhizosphere soil was defined as those shaken from the roots while rhizosphere soil was defined as those that adhere to the roots tightly and cannot be shaken off (Chaudhary et al., 2015). After a bunch of plants with soil was dug out, we shook the plants to collect the non-rhizosphere soil that departed from the plant roots until no soil could be shaken off. Then we brushed the roots to collect the rhizosphere soil until no soil could be brushed off (Chaudhary et al., 2015;Chen et al., 2020;Fu et al., 2022). The plants of C. songorica were separated into shoots and roots and oven-dried at 65°C to a constant weight. Soil samples from the same plot were mixed to form a composite sample; therefore, we obtained 12 rhizosphere and 12 non-rhizosphere composite soil samples (four treatments × thee replications).
Each composite soil sample was equally separated into three portions: one air-dried, one kept fresh at 4°C, and one stored at −80°C. The fresh sample was used to determine soil moisture and inorganic nitrogen concentration. Soil moisture was measured using the gravimetric method. The weighted fresh soil (approximately 20 g) was oven-dried at 105°C to a constant weight. Soil moisture was determined by the difference in weight between fresh soil and ovendried soil. Concentrations of ammonium nitrogen and nitrate nitrogen were measured using a continuous flow analyzer (AutoAnalyzer III; Bran+Luebbe), and the inorganic nitrogen was their sum. The airdried samples were used to determine the concentrations of total carbon, total nitrogen, total phosphorus, available phosphorus, and pH values. Concentrations of total carbon and total nitrogen in soil/ root samples were determined by an element analyzer (Elementar Vario Micro Cube). The air-dried soil samples were digested with a mixture of H 2 SO 4 -H 2 O 2 . Then the concentration of soil total phosphorus was determined by the molybdenum blue-ascorbic acid method on a spectrophotometer. The available phosphorus was determined by the same methods after extraction with 0.5 M NaHCO 3 .
Soil pH was measured using a pH meter. The frozen samples were used to extract DNA for evaluating the richness and composition of microbial communities.

| DNA extraction, PCR amplification, and sequencing
For each sample, DNA was extracted from 0.25 g of thawed soil using a DNeasy® PowerSoil® DNA isolation kit (Qiagen). The quality of isolated DNA was determined by electrophoresis on 1% agarose gel.
The purity was examined using NanoDrop. In total, 12 rhizosphere and 12 non-rhizosphere DNA samples were prepared.
The 16S rRNA gene was amplified to estimate the relative abundance of soil bacteria while internal transcribed spacer (ITS) region for the fungi. The paired primers, 341F and 806R, were used to target the hypervariable V3 + V4 region of the 16S rRNA gene; while the paired primers, ITS5-1737F and ITS2-2043R, were used for the ITS1 regions of fungi. The primers in each sample were tagged with unique barcodes. All PCRs were carried out in 30-μL reactions with 10 ng of template DNA, 0.2 μM of forward and reverse primers, and 15 μL of Phusion® High-Fidelity PCR Master Mix (New England Biolabs). For thermal cycles, the initial denaturation was at 98°C for 1 min, followed by 30 cycles of denaturation at 98°C for 10 s. Then annealing was at 55°C for 30 s and elongation at 72°C for 30 s. Finally, an extension was at 72°C for 5 min.
The PCR products of each sample and 1× loading buffer (with SYBR green) were mixed with equal volumes and then electrophoresed on a 2% agarose gel to examine the successful amplification of target DNA fragments. The PCR products for each sample were quantified and mixed with 80 ng of DNA in equal proportions.
The resulting mixed PCR products were electrophoresed with 2% agarose gel, and the DNA fragments ranging from 400 to 450 bp were excised and then purified with a GeneJET™ Gel Extraction Kit (Thermo Scientific).
Following the manufacturer's recommendation, the Ion Plus Fragment Library Kit 48 rxns (Thermo Scientific) was used for sequencing libraries preparation. The Qubit@ 2.0 Fluorometer (Thermo Scientific) and the Agilent Bioanalyzer 2100 system (Agilent) were used to evaluate the quality of the library. The libraries were sequenced on the Ion S5™XL platform, generating 400-bp/600-bp single-end reads.

| Bioinformatics analyses
Based on the unique barcode of each sample, single-end reads were identified and truncated by shearing off the barcode and primer sequences. The Cutadapt quality-controlled process was used to filter the raw reads under specific filtering conditions (Martin, 2011), we obtained the high-quality clean reads.
The reads were compared with the reference database (Silva and Unite database) using the UCHIME algorithm to detect chimera sequences (Edgar et al., 2011;Kõljalg et al., 2013;Quast et al., 2013) and remove them. Then, final clean reads were generated. We performed sequence analysis using Uparse software (Uparse v7.0.1001; Edgar, 2013). Sequences with a similarity >97% were assigned to the same operational taxonomic unit (OTU). The representative sequence for each OTU was screened for further annotation. The taxonomic information of bacterial representative sequences was annotated using the SILVA132 SSU rRNA database on Mothur algorithm (Quast et al., 2013). For fungal representative sequences, we annotated taxonomic information using the UNITE (v7.2) database-based BLAST algorithm (Kõljalg et al., 2013), which was calculated by QIIME software (Version 1.9.1).
In total, we obtained 947,987 and 958,264 bacterial high-quality sequences from rhizosphere and non-rhizosphere samples, respectively. These reads were clustered into 23,160 and 37,278 OTUs. For fungi, we obtained 961,649 and 928,290 high-quality sequences from rhizosphere and non-rhizosphere samples, respectively, which were clustered into 11,695 and 17,424 OTUs (Table A1).

| Statistical analysis
All statistical analyses were conducted in R 4.0.2 (R Core Team, 2020).
We used the OTU richness, Shannon diversity, ACE, Chao1, and phylogenetic diversity of bacterial and fungal communities to evaluate the response of microbial diversity. Before analyses, Shapiro-Wilk test was used for normal distribution evaluation, and Bartlett's test for variance homogeneity of the data. Differences in bacterial and fungal diversity between rhizosphere and non-rhizosphere zones through the Wilcoxon rank-sum test, performed using the wilcox.test function in the package 'stats' (Field et al., 2012). Difference in composition (e.g. the relative abundance of main bacterial and fungal phylum or class) between rhizosphere and non-rhizosphere zones was determined by the Wilcoxon rank-sum test. Furthermore, to identify significant differences in bacterial and fungal taxa between rhizosphere and nonrhizosphere zones, we used the linear discriminant analysis (LDA) effect size (LEfSe) method (Segata et al., 2011). The LEfSe analysis applied the Kruskal-Wallis sun-rank test to detect taxa with significantly different abundances between rhizosphere and non-rhizosphere zones and signed LDA log-score (LDA ≥ 4 in this study) to estimate the effect size of each different abundance taxon (Segata et al., 2011).
One-way analysis of variance (ANOVA) was used to assess the effects of grazing intensity on microbial diversity, composition, biotic variables (shoot biomass, root biomass, root total carbon, and root total nitrogen of C. songorica), and abiotic variables (concentrations of total carbon, total nitrogen, and total phosphorus, concentrations of available nitrogen and available phosphorus, soil pH, and soil moisture). After ANOVA, the Bonferroni test was used to perform a multiple comparison test that examined the differences in microbial diversity, composition, biotic and abiotic variables between treatments. Differences in microbial community composition between the two zones were evaluated using permutational multivariate analysis of variance (PERMANOVA) with the function "adonis" in the package "vegan" with 999 permutations (Oksanen et al., 2020), and visualized by Principal Coordinate Analysis (PCoA). The dissimilarity in microbial community composition was determined by the Bray-Curtis distance between treatments, performed with the "vegdist" function in package "vegan" (Oksanen et al., 2020).
A Mantel test was used to evaluate how changes in microbial community diversity or composition are related to biotic and abiotic variables by the "mantel test" function in the package "vegan", and visualized by the package "linkET" (Huang, 2021). The multicollinearity of the variables was assessed based on variance inflation factors (VIF), and those with a VIF value <10 were regarded as variables with low collinearity. The relative importance of variables in affecting microbial richness and composition was evaluated with the package of "rdacca.hp" (Lai et al., 2021).

| Microbial richness and composition in rhizosphere versus non-rhizosphere
Microbial communities differed markedly between rhizosphere and non-rhizosphere soils in OTU richness ( Figure 1). The bacterial richness was significantly lower in rhizosphere than in non-rhizosphere ( Figure 1a), and so was the fungal richness ( Figure 1b).
The bacterial community composition differed greatly between the two zones ( Figure 1c). This was more evident in the relative abundances of three dominant bacterial phyla, Pseudomonadota, Acidobacteria, and Actinobacteria. The relative abundance of Pseudomonadota was 74.4% in rhizosphere but decreased to 23.6% in non-rhizosphere ( Figure 1e; Table A2). In contrast, the relative abundances of Acidobacteria and Actinobacteria were 6.4% and 5.9% in rhizosphere but increased to 21.2% and 26.7% in non-rhizosphere, respectively ( Figure 1e; Table A2). Similarly, at the class level, the relative abundance of Gammaproteobacteria was 6.6% in non-rhizosphere, while it increased to 59.5% in rhizosphere ( Figure 1g; Table A3). By contrast, the relative abundances of Thermoleophilia, Acidobacteria, Acidimicrobiia, Actinobacteria, Deltaproteobacteria, Gemmatimonadetes, and Rubrobacteria were significantly higher in non-rhizosphere soil than in rhizosphere soil ( Figure 1g; Table A3).

F I G U R E 1
Comparisons in microbial community richness and composition between rhizosphere (Rs) and non-rhizosphere (Non-Rs) soils. Note: Bacterial (a) and fungal (b) richness; Difference in bacterial (c) and fungal (d) community composition; Relative abundance of main bacterial (e) and fungal (f) taxa at phylum level; Relative abundance of main bacterial (g) and fungal (h) taxa at class level. The taxa with a relative abundance <1% were assigned into "Others." *, **, and *** indicated significant difference between rhizosphere and non-rhizosphere soils by t-test at p = .05, .01, and .001, respectively.
The fungal community composition also differed between the two zones ( Figure 1d). For example, the relative abundance of Ascomycota was 48.8% in rhizosphere soil but decreased to 38.6% in non-rhizosphere soil (Figure 1f; Table A2). At the class level, the relative abundances of Dothideomycetes and Agaricomycetes were 24.6% and 12.2% in rhizosphere soil, while decreased to 10.41% and 4.01% in non-rhizosphere soil, respectively ( Figure 1h; Table A3).

| Effects of grazing intensity on microbial richness and composition
Grazing intensity had significant impacts on rhizosphere bacterial richness but not on that in non-rhizosphere soil ( Figure 2a). However, grazing intensity had little impact on fungal richness in both zones ( Figure 2b).
For the microbial composition, the dissimilarity in the bacterial community exhibited an increasing trend with the increase in grazing intensity in both zones (Figure 3a). In rhizosphere soil, the relative abundance of Pseudomonadota increased significantly under moderate and heavy grazing while that of Acidobacteria decreased ( Figure 3b; Table A4). In non-rhizosphere soil, the relative abundances of Pseudomonadota and Acidobacteria decreased with heavy grazing, while that of Actinobacteria increased (Figure 3c; Table A4). The dissimilarity in the fungal community also displayed an increasing trend with the increase in grazing intensity in both zones ( Figure 3d). Specifically, the dominant taxa, Ascomycota, showed a marked increase in relative abundance under moderate and heavy grazing (Figure 3e,f; Table A4).

| Biotic and abiotic drivers for microbial richness
In general, all the biotic and abiotic variables examined, except soil total phosphorus and soil pH, significantly affected soil bacterial richness in the rhizosphere, whereas that in the non-rhizosphere was not affected by these variables (Figure 4). The fungal richness in both zones was not affected by examined biotic and abiotic variables ( Figure 4).

| Biotic and abiotic drivers for microbial community composition
In general, all the biotic and abiotic variables examined, except soil total phosphorus and soil pH, affected the composition of the bacterial or fungal communities in both zones (Figure 5a). Moreover, among six variables with low collinearity, the rhizosphere bacterial composition was mainly affected by root total nitrogen and soil total nitrogen, whereas that in non-rhizosphere was mainly affected by soil-available phosphorus (Figure 5b,c). In addition, soil-available phosphorus was consistently the primary variable affecting the fungal composition in both rhizosphere and non-rhizosphere (Figure 5d,e).

| Microbial diversity and composition in rhizosphere and non-rhizosphere soils as influenced by grazing intensity
Our results demonstrated that the microbial communities differed markedly between rhizosphere and non-rhizosphere in diversity and composition in this desert steppe ecosystem. The microbial diversity indices in rhizosphere, including richness, Shannon diversity, ACE, Chao1, and phylogenetic diversity, were lower than those in non-rhizosphere (Figures 1 and A1), which is in line with previous findings that microbial diversity decreased with decreasing distance to plant roots (Donn et al., 2015;Fan et al., 2017), likely because the rhizosphere microbiomes are derived from the non-rhizosphere microbiomes under the selection effect of host plants (Hartmann et al., 2009;Nan et al., 2020).
The microbial composition also differed greatly between the two zones (Figure 1c,d). This was especially evident for the relative abundances of dominant microbial taxa. The relative abundance of Pseudomonadota, a dominant bacterial phylum, was 23.6% in non-rhizosphere but increased to 74.4% in rhizosphere (Figure 1e; Table A2). Further LEfSe analysis indicated that the increased abundance of Pseudomonadota can be attributed to the enrichments in Gammaproteobacteria, Pseudomonadales, Pseudomonadaceae, Pseudomonas, Pseudomonas frederiksbergensis at class, order, family, F I G U R E 2 Effects of grazing intensity on microbial community richness in rhizosphere (Rs) and non-rhizosphere (Non-Rs) soils. Note: G0, G0.15, G0.30, and G0.45 refer to control, light, moderate, and heavy grazing, respectively. Different lowercase letters indicate significant difference between treatments at p = .05. genus and species levels in rhizosphere zone under light or heavy grazing intensity ( Figure A4). Such changes are likely due to their similar strategy in nutrient use. Pseudomonadota is a representative copiotrophic bacteria phylum (Nie et al., 2018;Zhang et al., 2018. The increased availability of nutrients in rhizosphere due to acid activation by root exudates fosters its growth (Cesco et al., 2010), because more nutrients facilitate copiotrophic bacteria but restrict the oligotrophic bacteria according to the oligotrophcopiotroph hypothesis (Fierer et al., 2007;Leff et al., 2015). Similarly, Ascomycota is a dominant fungal phylum, and its relative abundance increased from 38.6% in non-rhizosphere to 48.8% in rhizosphere ( Figure 1f; Table A2). A similar result had been observed in an alpine grassland . Our further analyses indicated that the higher abundance of Ascomycota in rhizosphere zone mainly stemmed from the higher abundances of Dothideomycetes and Agaricomycetes at the class level (Figure 1h). The relative abundances of Dothideomycetes and Agaricomycetes increased from 10.41% and 4.01% in the non-rhizosphere soil to 24.55% and 12.16% in the rhizosphere soil, respectively (Table A3). Collectively, these results suggest that the microenvironment of the plant rhizosphere may facilitate the flourishing of certain bacterial and fungal taxa.
An intriguing result from our experiment was that the rhizosphere bacterial diversity was more sensitive to grazing intensity than that in non-rhizosphere . As indicated in Figure 2a, the rhizosphere bacterial richness declined dramatically under heavy grazing intensity, while that in non-rhizosphere remained unchanged . Such differing responses also held for other diversity indices, including Shannon diversity, ACE, Chao1, and phylogenetic diversity ( Figure A2). These consistent results suggest that grazing intensity has a strong impact on rhizosphere bacterial diversity. Two mechanisms may underlie the grazing-induced decline in rhizosphere bacterial diversity. First, herbivory removes the aboveground organs of a plant, resulting in a reduction in the assimilated carbon allocated to belowground organs, and thus a decline in the quantity of carbon sources available to bacteria in the rhizosphere of the plant (Aldezabal et al., 2015;Byrnes et al., 2018;Mueller et al., 2017), which in turn reduces rhizosphere bacterial diversity. In this study, the rhizosphere bacterial richness was significantly affected by root carbon (Figure 4), providing support for this mechanism. Second, plant diversity is usually coupled with microbial diversity in a community (De Deyn & van der Putten, 2005). Previous studies have suggested that biodiversity F I G U R E 3 Effects of grazing intensity on the dissimilarity of microbial community composition and on the relative abundance of microbial taxa in rhizosphere (Rs) and non-rhizosphere (Non-Rs) soils. Note: The dissimilarity of bacterial communities (a) and the dissimilarity of fungal communities (d); The relative abundance of bacterial phyla in rhizosphere (b) and in non-rhizosphere (c); Relative abundance of fungal phyla in rhizosphere (e) and in non-rhizosphere (f). Different capital and lowercase letters indicate significant differences between treatments in non-rhizosphere and rhizosphere soils, respectively. G0, G0.15, G0.30, and G0.45 are the same as in Figure 2. can be enhanced by intermediate disturbance (Connell, 1978).
However, in this desert steppe, long-term grazing, even under low intensity, greatly reduced the species diversity .
Grazing-induced reduction in plant diversity can decrease the diversity of carbon sources for bacterial taxa (Dwivedi et al., 2017), as a consequence, some bacterial taxa cannot grow well or disappear and the rhizosphere bacterial diversity declines.
In contrast to the response of rhizosphere bacterial richness, the fungal diversity was not significantly affected by grazing intensity in both zones (Figures 2b and A2), which is consistent with results from grazing experiments Zhang & Fu, 2020). The weak response of fungal richness may be due to the fact that fungi were more tolerant to environmental stresses than bacteria (Rousk & Bååth, 2011). For example, fungi often exhibit more resistance to soil drought or acidity than bacteria (Rousk & Bååth, 2011). In this study, the fungal richness was not significantly affected by the biotic and abiotic factors examined (Figure 4), implying that the grazinginduced changes in these variables did not exert a significant impact on the fungal richness. Overall, our results suggest that grazing on this grassland has little impact on the soil fungal richness.
Our results suggest that grazing intensity has stronger impacts on bacterial composition in rhizosphere than that in non-rhizosphere.
For example, heavy grazing resulted in a 6.55-fold increase in community dissimilarity in rhizosphere while a 1.89-fold increase in non-rhizosphere (Figure 3a). An impressive result of this experiment was that, with the increase of grazing intensity, the relative abundance of a dominant bacterial phylum, Pseudomonadota, decreased in non-rhizosphere but increased in rhizosphere (Figure 3b,c). Such opposite responses to grazing intensity were likely due to its copiotrophic strategy in nutrient use (Leff et al., 2015). The increased grazing intensity led to a relatively higher plant root nitrogen (Tables A6 and A7), which in turn increased the relative abundance of Pseudomonadota in rhizosphere due to the positive relationship between the two ( Figure A5a). By contrast, the grazing-induced decrease in relative abundance of Pseudomonadota in non-rhizosphere was related to the reduction in soil total nitrogen ( Figure A5b; Tables A6   and A7). This was especially the case for Alphaproteobacteria, as its relative abundance was closely related to soil nitrogen ( Figure A6a; Tables A6 and A7), and a significant decrease in Alphaproteobacteria was observed in non-rhizosphere zone under moderate and heavy grazing intensities ( Figure A3b; Table A5).
Our results also demonstrated strong impacts of grazing intensity on fungal composition in this grassland. As indicated by the increased dissimilarity of the fungal community with grazing intensities (Figure 3d), the fungal composition was dramatically altered by heavy grazing in both zones. A significant increase in relative abundance was observed for a dominant taxon, Ascomycota, under heavy grazing intensity (Figure 3e,f). This phenomenon could be explained by its oligotrophic strategy, which facilitates its growth in a relatively low-carbon environment (Sterkenburg et al., 2015). As indicated in this experiment, heavy grazing resulted in a great reduction in carbon concentrations in plant roots and soils (Tables A6   and A7), and the relative abundance of Ascomycota was negatively related to these carbon concentrations ( Figure A5c,d).

F I G U R E 4
Impacts of abiotic and biotic factors on bacterial and fungal richness in rhizosphere (Rs) and non-rhizosphere (Non-Rs) soils. RBI, root biomass; RTC, root total carbon; RTN, root total nitrogen; SAN, soil-available nitrogen; SAP, soil-available phosphorus; SBI, shoot biomass; SM, soil moisture; STC, soil total carbon; STN, soil total nitrogen; STP, soil total phosphorus. Overall R 2 indicates the total explained fraction by the six variables.

| The roles of biotic and abiotic factors in affecting microbial responses to grazing intensity
In this experiment, the grazing-induced decrease in rhizosphere bacterial richness was related to changes in 9 of 11 variables ( Figure 4), consistent with a previous study in alpine meadow finding that soil microbial richness was affected by multiple biotic and abiotic factors (Zhang & Fu, 2020). Moreover, we found that all the examined biotic variables were closely related to rhizosphere bacterial richness (Figure 4), highlighting the impacts of plant selection F I G U R E 5 Impacts of abiotic and biotic factors on microbial compositions in rhizosphere (Rs) and non-rhizosphere (Non-Rs) soils. Note: Mantel results on the impacts of the variables examined on bacterial and fungal compositions (a). Relative roles of six variables with low collinearity on rhizosphere (b) and non-rhizosphere bacterial compositions (c). Relative roles of six variables with low collinearity on rhizosphere (d) and non-rhizosphere fungal compositions (e). SBI, RBI, RTC, RTN, STC, STN, STP, SAN, SAP, and SM are the same as in Figure 4. Overall R 2 indicates the total explained fraction by the six variables. effect on rhizosphere bacterial richness (Hartmann et al., 2009;Nan et al., 2020). However, the examined variables had little effect on fungal richness in both zones (Figure 4), consistent with previous results (Yang, Niu, et al., 2018;Yang, Zhu, et al., 2018). This may be due to the fact that fungi are more resistant to environmental variability relative to bacteria (Rousk & Bååth, 2011).
Previous studies have highlighted the role of soil pH in affecting microbial richness (Qu et al., 2016;Yang, Niu, et al., 2018;Yang, Zhu, et al., 2018;Zhang & Fu, 2020). However, in this study, the changes in bacterial and fungal richness were not affected by soil pH in both zones. The weak effect of soil pH in our study could be explained by the results that soil pH in this experiment was not significantly altered by grazing intensity (Tables A6 and A7).
Our results indicated that the bacterial and fungal compositions in both rhizosphere and non-rhizosphere zones were co-affected by several biotic and abiotic variables (Figure 5a). The rhizosphere bacterial composition was primarily affected by the root biomass and the root carbon concentration, highlighting the importance of plant-related resources in affecting rhizosphere bacterial communities (Berg & Smalla, 2009;Costa et al., 2006). In contrast, the bacterial composition in non-rhizosphere soil was primarily affected by soil-available phosphorus, highlighting the role of abiotic factors in affecting non-rhizosphere bacterial composition (Schöps et al., 2018;Tian et al., 2017). Moreover, soil-available phosphorus was also the main factor that affected the fungal composition in both zones (Figure 5d,e). The important role of soil available-phosphorus in influencing the microbial composition may be due to the fact that grasslands in northern China have relatively lower soil phosphorus than other grasslands (Han et al., 2005). For example, the soil total phosphorus concentration (515 mg kg −1 ) in this experiment site (Table A7) was much lower than the mean value (699 mg kg −1 ) in grasslands of North America (US Geological Survey, 2001). In addition, soils in this region are rich in Ca 2+ and Mg 2+ , and phosphorus is usually bounded to these metal ions, further constraining the available phosphorus for soil fungi . Importantly, the lack of soilavailable phosphorus was exacerbated by heavy grazing in this experiment (Tables A6 and A7). This is likely due to grazing-induced losses of fertile soils. The Mongolia steppes are located in the west wind belt and are characterized by strong winds throughout the year.
Years of heavy grazing caused a decrease in vegetation cover and an increase in bare land patches. As a result, strong winds caused severe wind erosion, such as sand-dust storms, by which the nutrient-rich soils in top layers were removed from grasslands (Giese et al., 2013).
Given that some arbuscular mycorrhizal fungi can help plants efficiently exploit and absorb soil phosphorus (Faghihina et al., 2020;van der Heijden et al., 1998;van der Heyde et al., 2017), the strong relationship between soil-available phosphorus and the fungal community composition also suggests that some fungal taxa may play an important role in mediating the availability of soil phosphorus in this grassland. Overall, our study provides experimental evidence that grazing intensity caused changes in abiotic and biotic factors, and these changes in turn altered the soil microbial composition.

CO N FLI C T O F I NTER E S T S TATEM ENT
The authors declare that they have no competing interests.

O PE N R E S E A RCH BA D G E S
This article has earned an Open Data badge for making publicly available the digitally-shareable data necessary to reproduce the reported results. The data is available at https://github.com/IBCAS yangy ang/Publi c-data-for-ECE.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are available from Note: Main bacterial and fungal phyla indicate the phyla with relative abundance higher than 1%.

TA B L E A 3
The relative abundances of main bacterial and fungal class under different grazing intensities in rhizosphere and non-rhizosphere zones. Note: Significant effects are highlighted in bold (p < .05). The soil zones refer to plant rhizosphere and non-rhizosphere.

TA B L E A 6
Results (F and p values) of one-way ANOVA of the effects of grazing intensity on biotic and abiotic variables.

Variables
Grazing intensity treatment F I G U R E A 2 Effects of grazing intensity on soil bacterial and fungal diversity indices in rhizosphere and non-rhizosphere soils. Note: G0, G0.15, G0.30, and G0.45 refer to control, light, moderate, and heavy grazing, respectively. Different lowercase letters indicate significant difference among treatments at p = .05.

F I G U R E A 3
Effects of grazing intensity on the relative abundance of main soil bacterial and fungal classes in rhizosphere (a, c) and non-rhizosphere (b, d) soils. Note: G0, G0.15, G0.30, and G0.45 refer to control, light, moderate, and heavy grazing, respectively. Different uppercase and lowercase letters indicate significant difference among grazing intensities in rhizosphere and non-rhizosphere soils, respectively. The abundance of a letter indicates that the difference is not significant among grazing intensities in both soils. Significant level is p = .05.

F I G U R E A 4
Identity biomarkers with linear discriminant analysis (LDA) scores of 4 or greater in bacterial community between rhizosphere (Rs) and non-rhizosphere (Non-Rs) soils under light (a) and heavy (b) grazing intensities. Note: LEfSe analysis (LDA ≥ 4) only detected identify biomarkers between rhizosphere and non-rhizosphere zones under light and heavy grazing intensities, while no identify biomarkers of fungi were detected between two zones under all grazing intensities.

F I G U R E A 5
Relationship of the relative abundance of main bacterial phylum in non-rhizosphere (a) and in rhizosphere (b) soil, main fungal phylum in non-rhizosphere (c) and in rhizosphere (d) soil with biotic and abiotic variables. ***, p < .001; **, p < .01; *, p < .05. RBI, root biomass; RTC, root total carbon; RTN, root total nitrogen; SAN, soil-available nitrogen; SAP, soil-available phosphorus; SBI, shoot biomass; SM, soil moisture; STC, soil total carbon; STN, soil total nitrogen; STP, soil total phosphorus.

F I G U R E A 6
Relationship of the relative abundance of main bacterial class in non-rhizosphere (a) and in rhizosphere (b) soil, main fungal class in non-rhizosphere (c) and in rhizosphere (d) soil with biotic and abiotic variables. ***, p < .001; **, p < .01; *, p < .05. RBI, root biomass; RTC, root total carbon; RTN, root total nitrogen; SAN, soil-available nitrogen; SAP, soil-available phosphorus; SBI, shoot biomass; SM, soil moisture; STC, soil total carbon; STN, soil total nitrogen; STP, soil total phosphorus.